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The grand potential for open systems describes thermodynamics of fluid flows at low Mach num- 
bers. A new system of reduced equations for the grand potential and the fluid momentum is derived 
from the compressible Navier-Stokes equations. The incompressible Navier-Stokes equations are 
the quasi-stationary solution to the new system. It is argued that the grand canonical ensemble is 
the unifying concept for the derivation of models and numerical methods for incompressible fluids, 
illustrated here with a simulation of a minimal Boltzmann model in a microflow setup. 
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The classical incompressible Navier-Stokes equation 
(INS) is mechanical description of fluid flows at low Mach 
numbers (Mach number, Ma = Uo/c s , is the ratio of the 
characteristic flow velocity Uo to the isentropic sound 
speed c s denned at some reference temperature T and 
density p). The INS equation can be written in the Eu- 
lerian coordinate system as, 

d t u a + UpdpUa + d a P = — dpdpu a , d a u a = 0, (1) 

where u is the fluid velocity, P is the pressure and Re 
the Reynolds number, which characterizes the relative 
strength of the viscous and the inertial forces Q. The 
pressure in Q is not an independent thermodynamic 
variable but is rather determined by the condition of the 
incompressibility, 

dpdpP = - (dpu a ) (d a up) . (2) 

Thus, in order to obtain the pressure at a point, one has 
to solve the Laplace equation (0) in a domain, and the re- 
lationship between the pressure and the velocity becomes 
highly nonlocal. The physical meaning of J2J) is that in a 
system with infinitely fast sound propagation, any pres- 
sure (and thus density) disturbance induced by the flow 
is instantaneously propagated into the whole domain. 

A textbook justification for the thermodynamics of the 
INS description is usually based on the isentropic flow as- 
sumption. If this assumption is valid for all times (the 
entropy density is simply convected by the flow), then the 
thermodynamic pressure depends only on the acoustic 
(density) variations. For low Mach number flows, these 
variations adjust to the flow on every spatial scale in 



the long time dynamics [lj,|2|. The way this adjustment 
takes place is nontrivial, and it was given considerable 
attention recently jj, U, L| In particular, it was proved 
that weak solutions of the isentropic compressible Navier- 
Stokes equations converge to that of the INS equation for 
some special boundary conditions (adiabatic absorbing 
walls) 5]. However, it remains a challenge to give a ther- 
modynamic derivation of the incompressibility without 
the isentropic flow assumption. Such a thermodynamic 
derivation goes far beyond academic interest. Indeed, as 
is well known, the INS equations are extremely hard to 
study, both analytically and numerically. Therefore, an 
extended system where the flow is coupled to a dynamic 
equation for a scalar thermodynamic variable provides a 
better starting point for numerical and theoretical stud- 
ies of the incompressible hydrodynamics. Indeed, in the 
computational fluid dynamics, at least two undeniably 
successful routes to avoid the "elliptic solver problem" , 
that is, avoiding the nonlocality of the pressure @ are 
well known. The first is the so-called artificial compress- 
ibility method introduced by Chorin [(| and Temam 0, 
where an evolution equation for the pressure is postu- 
lated instead of the constraint (JSJ (see e. g. |§| for a 
recent review). The second route is kinetic-theory mod- 
els such as the lattice Boltzmann method 0, where a 
Boltzmann-like equation is obtained for the compressible 
fluid flow in the low Mach number limit. The thermo- 
dynamics of the lattice Boltzmann method was clarified 
|lOj, and the method enjoys a thermodynamically sound 
derivation from the Boltzmann equation Can the 

compressibility methods be modified in a way to make 
them physical models? Furthermore, for emerging fields 
of fluid dynamics such as flows at a micrometer scale, 
corrections to the incompressibility assumption become 
crucial [l2j |. 
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In this work, we present the thermodynamic descrip- 
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tion of incompressible fluid flows. In step one, we will 
argue that the grand potential is the proper thermo- 
dynamic potential to study the onset of incompressibil- 
ity. Use of the grand potential instead of the entropy 
enormously simplifies the equations of compressible hy- 
drodynamics in the low Mach number limit. We show 
that after a short time dynamics during which acous- 
tic waves are damped by viscosity, the fast dynamics of 
the grand potential becomes singularly coupled to the 
slow dynamics of momentum, and reduced compressible 
Navier-Stokes equations are derived (RCNS). The incom- 
pressible Navier-Stokes equation is the quasi-stationary 
solution of the RCNS, when Mach number tends to 
zero. RCNS equations are a thermodynamically consis- 
tent generalization of the compressibility schemes. Fi- 
nally, by writing the grand potential for the Boltzmann 
equation, we show that, upon an appropriate discretiza- 
tion of velocities, the present constructions leads to the 
cntropic lattice Boltzmann method. Correctness of the 
present almost-incompressible description is illustrated 
with a simulation of a microflow setup. 

A low Mach number flow is a setup where only small 
spatial deviations of the entropy and the density from 
the equilibrium value exist. The grand potential is the 
natural thermodynamic variable to describe such a setup 
(and the corresponding to it grand canonical ensemble is 
the natural statistical thermodynamics framework). This 
happens because the balance laws (compressible Navier- 
Stokes equations or the Boltzmann equation) are always 
written for a sufficiently small volume element in the 
Eulcrian frame of reference (i.e. the volume element is 
fixed in space). From a thermodynamic standpoint, this 
volume element is an open system. From elementary 
thermodynamics we know that in an open system, the 
thermodynamic equilibrium is conveniently described in 
terms of the grand potential Q(tp,T), where ip is chemical 
potential, and T is the temperature. 



Now, we shall find the expression for the grand po- 
tential in the Eulerian coordinate system for a vol- 
ume element SV, in thermodynamic equilibrium. In 
the co-moving system, the grand potential is written as 
n L (V>,T) = -P{ip,T)5V, where P is the pressure. The 
transition to the Eulerian (fixed) system is done by fixing 
the momentum, Oe(V') T, m) = — PSV + X a m a SV, where 
X a are Lagrange multipliers, and m is the momentum 
density. For small values of momentum, the Lagrange 
multipliers can be specified by noting that the energy in 
the Eulerian coordinate system is [e+(m 2 /2p)]SV, where 
e is the internal energy density. Using the relationship be- 
tween the energy and the grand potential for the thermo- 
dynamic equilibrium, we find that X a = m a /2p + 0(m 3 ). 
Thus, in the Eulerian coordinate system, the grand po- 
tential, up to the higher-order terms in momentum, is 
written as, 



E = (-P + — j SV. (3) 

The difference of the pressure and the kinetic energy is 
the (negative of) density of the grand potential, and it 
will be used below as the natural thermodynamic poten- 
tial for the low Mach number flows: 



Dynamic equations for the set of variables p, m and 
Q are written using the standard compressible Navier- 
Stokes equations (CNS) for Newtonian fluids Note 
that CNS are usually written in terms of a different set 
of variables (for example, in terms of the entropy density 
S instead of Q). The recomputation from either form of 
the CNS equations to the present set of variables poses 
no difficulties, and we here write the final result: 



d t m a + dp 



G + 2L) Sap 



dtp + d a m a = 0, 
m a m,3 t rr = (i 



dP 
dp 




U a0 d a ?») 



(5) 



1 dP 



p ) pCydT 



d a (Kd a T) , 



where LI Q ^ is the stress tensor of a Newtonian fluid, 



n Q 3 = -n 



9/3 { — ) - <W f - a ] a 



with D the spatial dimension, p the shear viscosity, and 
A the ratio of bulk to shear viscosity. In (JSJ, k is the 



ermal conductivity, Cy is the specific heat at constant 
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volume, and the temperature T is known from the equa- 
tion of state. Variations of material parameters such as 
viscosity p will be ignored in the further discussion. Note 
that equations iJSJ are "exact" in the sense that they are 
just the standard CNS equations written for p, m, and 
the function Q However, the physical meaning of the 
function Q J3J as the density of the grand potential is 
valid only up to the lowest order in momentum, and thus 
the form of the CNS (0) is an important intermediate step 
in the study of the low Mach number hydrodynamics. 

Since the speed of sound, c s = \JdPjdp\s, contributes 
to the dynamic equation for the grand potential, it is 
instructive to rewrite © as: 



d t G + d a (c 2 m a ) -d a {pc 2 s ) 

P 



(6) 



We expect the following scenario of the onset of incom- 
pressibility, as it can be inferred from ©: If the speed of 
sound is "large" , then, after a "short-time" dynamics of 
the density leading tops; const, the third term in © can 
be neglected, whereas the second term becomes the dom- 
inant contribution to the time derivative of Q. The dy- 
namic equation for the grand potential becomes then sin- 
gularly perturbed, and represents the "fast" mode cou- 
pled to the "slow" dynamics of momentum. The contri- 
butions to the time derivative of the grand potential not 
displayed in © are responsible for corrections to the in- 
compressibility. We now proceed with quantifying these 
statements. 

For small deviations from the no-flow situation, 
(|m Q m/5| < P/i, which implies Ma <C 1), we are in- 
terested in the long time solutions of the CNS equations 
(times of the order of the momentum diffusivity time 
tmd ~ P~L 2 / [i, where L is a characteristic length associ- 
ated with the flow). We define a dimensionless number, 
Kn = p/(pc s L), the ratio of the sound propagation time 
L/cs and the momentum diffusivity time, as the Knud- 
sen number for a general fluid, and we are considering 
Kn -c Ma <C 1. For the sake of simplicity, we assume 
that the Prandtl number Pr ~ 1 in the subsequent anal- 
ysis. The short time dynamics is isentropic and linear, 
and it is well known that, away from boundaries, any 
density perturbation at a distance r away from the dis- 
turbance source decays as (see Ij, p. 300): 



dp (r, t) oc (LarL) 



-1/2 



exp 



(r - c s tf 
2La rL 



(7) 



Here a new dimensionless number La is defined as: 



L, = K»|2-* + a) + *^, (8) 



where 7 the ratio of specific heat at constant pressure and 
volume. La generalizes the notion of Knudsen number 
for an arbitrary fluid (we call it the Landau number in 
the honor of Landau, who explained the relevance of this 
number in the context of acoustic damping 0, p. 300). 
Thus, the short term hydrodynamics reveals the following 
length scale L a and the time scale t a (since we assume 
Pr ~ 1, we need not distinguish between La and Kn for 
the present purpose): 



7 Kn 



L 



7 KnL, 



(9) 



At the time scale larger than i a , and on the spatial scale 
larger than L a , the density of the fluid can be safely 
treated as a constant (in the usual isentropic theory, the 
characteristic time for the onset of incompressibility is of 
the order L/c s ^> t a ). Note that the length scale L a was 
also found in the derivation of the sub-grid model from 
kinetic theory 13 . On the time-space scale larger than 
(|9"|). we can neglect the density variation, and the tem- 
perature variation 5T (from the globally uniform value 
T) becomes a function of the grand canonical potential, 



ST: 



dT 
dP 



2p 



Once the time and space scales (0 are identified, we 
complete the reduction of the CNS equations © by 
merely rescaling the variables. The momentum is scaled 
by the characteristic momentum ~p~Uo (known from the 
initial or boundary condition), j = m/(p[/o), and we 
introduce the reduced grand canonical density, O = 
G/(j>Uq). Making time dimensionless with t a , space with 
L a , (t — > t/t a , x — > x/i a ), neglecting variations of the 
density, and taking into account the thermodynamic re- 
lation for the temperature mentioned above, the two last 
equations in Q reduce to the following scale- independent 
closed set of equations for the dimensionless grand poten- 
tial and momentum: 
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d t j a =-Madp 



* 9 = "Ma" 9 ' 



jajp + S + 



Ma 2 j 2 



ja 1 



+ VKn [ 1 + A - — ) dadfjjfs + VKndpdpjc 



(10) 



dP 
df 



— jijT ( 9 "J/3 + ^Ja) uA,./. + S/jja) - U + A - ) ; (./ ,/;„./„ i 



Note that all "material parameters" appearing in llOjl 
(CV, 7, A, k) are evaluated at equilibrium at p and T. 

The reduced set of compressible Navier-Stokes equa- 
tions (|1C)|> is valid for Kn <C Ma <C 1, and, as we ex- 
plained above, on the scales larger than acoustic scales 
©. Roughly speaking, (|10fl is what survives from the 
compressible Navier-Stokes equations just before the in- 
compressibility sets on. Indeed, the time derivative of 
becomes singularly perturbed as the Mach number tends 
to zero, and we recover the incompressibility condition, 
d a ja = 0, as the quasi-stationary solution of the system 
(|10fl . This solution, when substituted in the momentum 
equation, recovers the INS equation JQl with the usual 
accuracy of the order 0(Ma 2 ). Note that the velocity in 
the INS equations recovered from I jlOjl is u = j/p. Cor- 
rections to the quasi-stationary solution can be found in a 
systematic way (see, e. g. an d we do not address 

this here. The following point needs to be stressed: The 
dissipation terms (proportional to VKn) cannot be ne- 
glected in the equation for the grand potential I jlOp and 
simultaneously kept in the momentum equation. This 
is at variance with the artificial compressibility method 
0,0- In other words, the RCNS is the minimal thermo- 
dynamic system for incompressible hydrodynamics. 

In this Letter, we reported a new basic physical fact of 
fluid dynamics: Grand potential for low Mach number 
flows gives the thermodynamic description of the incom- 
pressible phenomena. The resulting system H10JI includes 
a local (non-advected) equation for the scalar thermody- 
namic field. What follows from this fact? Let us list 
some of the consequences: 

• With the corrections mentioned above, the numer- 
ical schemes of the artificial compressibility kind 
become a firm status of physical models. 

• The structure of the coupling between the flow and 
the thermodynamic variable hints at that the true 
incompressible flows are attractors of system (|10|l . 
with the INS as the leading-order approximation. 
Because of the singular perturbation nature, it may 
be easier to study attractors of (|10|) , rather than of 
the INS equations. 



System <|10|) can be a starting point for a systematic 
derivation of nonlinear models for heat transport 
in the nearly-incompressible fluids such as multi- 
phase fluids, polymeric liquids and melts etc (this 
is relevant even in the linear case, see, e. g. 



cjthi 
0). 



• In the celebrated Kolmogorov theory, equation for 
the kinetic energy is used to make predictions about 
the structure of the fully developed turbulence. 
The present thermodynamic approach unambigu- 
ously delivers the density of the grand potential as 
the scalar field associated with the incompressible 
fluid flow, and thus can be relevant to develop the- 
ories of turbulence through studying the resulting 
balance equation. 

Let us dwell on the use of (|10|) for numerical simu- 
lation of incompressible flows. Recall that the spectral 
methods [13 for the INS JTJ based 

on the Fast Fourier 
Transform (FFT) for solving the Laplace equation @ 
are very efficient for high Reynolds number flows in sim- 
ple geometries. On the other hand, thanks to a rela- 
tively simple structure of the system (|10[1 (the equation 
for the scalar variable contains no convected derivative), 
it can be addressed by a host of discretization methods 
(see e. g. ||). The system (I1U|I can be useful for simula- 
tion of flows in complex geometries and especially non- 
stationary problems, that is, in the situations where one 
seeks to avoid solving the Laplace equation by restoring 
to some relaxation schemes (artificial compressibility, lat- 
tice Boltzmann etc). Special attention should be payed to 
the fact that the system (|10|) contains terms of different 
order of magnitude. This situation is typical for all relax- 
ational schemes such as, for example, the lattice Boltz- 
mann method, and we do not discuss here how to deal 
with this issue numerically (see, e. g. 0,0, El)- Note, 
however, the important smoothing effect of the diffusion 
term (inversely proportional to the Prandtl number) in 
the equation for 9 which makes the physical system (fTTTf) 
more amenable to numerics, that is, less numerically stiff 
than artificial compressibility methods. These numerical 
issues remain out of the scope of the present Letter, and 
will be addressed in a separate publication. 

While the new system of reduced compressible Navier- 
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the incompressibility per se, it indicates that grand po- 
tential may be a relevant thermodynamic variable in flow- 
ing systems way beyond "ordinary" fluids. 
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FIG. 1: Velocity profile of the 2D body force driven Poiseuille 
flow at Kn = 0.035 and Ma = 0.01. 

Stokes equations (fTTTfl is our central result, we con- 
clude this Letter with a more general statement that 
the grand potential (and the relevant grand canoni- 
cal ensemble) can be implemented for eventually any 
more microscopic setup (and not only for compressible 
Navier-Stokes equations as above). This viewpoint pro- 
vides a unified setting for derivations of a variety of 
mesoscopic or molecular dynamics models for numerical 
simulation of incompressible and nearly-incompressible 
flows. As an illustration here, let us assume Boltzmann's 
description with the one-particle distribution function 
/(v,x, t). Starting from the general form of grand po- 
tential, Q(f,a,(3,X) = J f [in / + a + \ a v a + ftv 2 ] dv, it 
is easy to show that at equilibrium / eq (v, a, (3, A) (de- 
fined from Sg = 0), we have £(/ cq ) = £ cq (a,,3) + 
(A 2 /2) / / oq (v,a,/3,0)(w 2 /2)dv for small A. When the 
velocity integral £/(/, a, (3, A) is evaluated with the Gauss- 
Hermit quadrature with the weight exp(/3v 2 ) at fixed j3, 
one obtains the entropy function of the lattice Boltzmann 
method [Iollll| . The present alternative derivation based 
on the grand potential is new. Thus, the entropic lattice 
Boltzmann method is a valid physical model for nearly- 
incompressible flows, and can be used for finite but small 
Knudsen low Mach number flow problems often encoun- 
tered in the microflows I n Fig. \T\we present ex- 
cellent comparison between the analytical solution to the 
Bhatnagar-Gross-Krook (BGK) kinetic equation and 
the entropic lattice Boltzmann simulation for the body- 
force driven 2D Poiseuille microflow. A detailed study of 
microflows within the entropic lattice Boltzmann setting 
is presented elsewhere [lif . 

While this paper was in the revision process, we 
learned about a very recent paper |20j | which states the 
usefulness of the grand potential in the context of plasma 
turbulence simulations. Although [ijj does not deal with 
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